Purpose
Introduction
Setup
Begin by downloading and opening the following packages into your library. The gelnet library will be used to train the model, dplyr and gdata are used for data manipulation.
library(gelnet)
library(dplyr)
library(gdata)
Data
The the pcbc data is organized in a data frame with 99 samples as columns and 219 methylation probes as rows.
load("pcbc.data.Rda")
pcbc.data
load("pcbc.pd.f.Rda")
head(pcbc.pd.f)
Then find the mean center by subtracting the mean of each probe from the entire pcbc data. The mean of each probe just be in a numeric vector the same size as the number of probes, in this case 219.
m <- apply(pcbc.data, 1, mean)
m[1:5]
cg02927655 cg15948871 cg17676824 cg25552705 cg21434327
0.4130770 0.2652034 0.3821945 0.4168715 0.3776927
pcbc.data.2 <- pcbc.data - m
pcbc.data.2
Identify stem cells and break up all samples into 2 groups:
- Stem cell (‘X.tr’ object)
- not stem cell (‘X.bk’ object)
# Define PCBC groups (SC and non.SC)
M1_smp <- pcbc.pd.f[pcbc.pd.f$Diffname_short %in% "SC",] #SC
M2_smp <- pcbc.pd.f[!(pcbc.pd.f$Diffname_short %in% "SC"),] #non-SC
# Select PCBC data
X.tr <- pcbc.data.2[, as.character(M1_smp$UID)] # 44 samples
X.bk <- pcbc.data.2[, as.character(M2_smp$UID)] # 55 samples
Now we can begin to train the the one-class model with the gelnet function. The gelnet function can be used for Linear Regression, Binary Classification and One class Problems by using an iterative method called coordinated descent (Sokolov et al. 2016).
It has four main arguments described below:
- X: n by p matrix => transpose(‘X.r’)
- y: ‘NULL’ for one class models
- l1: coefficient for the L1-norm penalty => ‘0’
- l2: coefficient for the L2-norm penalty => ‘1’
Make sure you transpose the matrix so that the genes are listed as rows and samples as columns. Then store the results as a tsv file (pcbc-stemsig.p219.rda).
## Train a one-class model
mm <- gelnet(t(X.tr), NULL, 0, 1) #NULL for a one-class task
Training a one-class model
Iteration 1 : f = 0.6931472
Iteration 2 : f = 0.3320004
Iteration 3 : f = 0.3263175
## Store the signature to a file
save( mm, file = "pcbc-stemsig.p219.Rda")
Leave One Out Cross-Validation
To test how the model’s performance by using leave one out cross-validation. This process has three steps:
- Train model on non-left-out data
- Score the left-out sample against the background
- AUC = P( left-out sample is scored above the background )
# Cross-validation with linear model:
# Perform leave-one-out cross-validation
auc <- c()
for(i in 1:ncol(X.tr)) {
## Train a model on non-left-out data
X1 <- X.tr[,-i]
X1 <- as.matrix(X1)
K <- t(X1) %*% X1 / nrow(X1)
m1 <- gelnet.ker(K, NULL, lambda = 1)
w1 <- X1 %*% m1$v
## Score the left-out sample against the background
X.bk <- X.bk[rownames(X.tr),]
X.bk <- as.matrix(X.bk)
s.bk <- t(w1) %*% X.bk
s.bk <- unmatrix(s.bk)
s1 <- t(w1) %*% X.tr[,i]
s1 <- unmatrix(s1)
## AUC = P( left-out sample is scored above the background )
auc[i] <- sum(s1 > s.bk) / length(s.bk)
cat( "Current AUC: ", auc[i], "\n" )
cat( "Average AUC: ", mean(auc), "\n" )
}
If the validation is successful you will notice that the auc variable will be a numeric vector consisting of only 1’s.
r1:c1
4.730943
[1] 1 1 1 1 1 1
[1] TRUE
Replace NA
The Replace NA function is used to replace any values that are NA (not available) with either the mean or the median value of the probe for a give group.
- check for NA values
- locate the NA values
- calculate the mean or median for the probe where each NA is found
- replace the values
replace.NA <-function(data,type.info,by = "mean"){
if(!"group" %in% colnames(type.info)) stop("type.info must have group column")
if(!"sample" %in% colnames(type.info)) stop("type.info must have a sample column")
# Do we have NAs?
if(is.na(table(is.na(data))["TRUE"])){
message("No NAs were found")
return(data)
}
# get NAs index
idx <- which(is.na(data) == TRUE,arr.ind=TRUE)
count <- table(rownames(idx))
message("======= Status Number of NA in probes ========")
message("--------------------- Summary------------------")
print(summary(as.numeric(count)))
message("\n----------- Probes with more nb of NAs -----------")
print(head(sort(count,decreasing = T)))
message("===============================================")
idx <- cbind(idx, mean = NA, median = NA)
# For each NA value calculate the mean for the same probe for the samples
# where it belongs
for(line in 1:nrow(idx)){
row <- idx[line,1]
col <- idx[line,2]
probe <- rownames(idx)[line]
sample <- colnames(data)[col]
group <- type.info[type.info$sample == sample,"group"]
samples.in.group <- type.info[type.info$group == group,]$sample
# get the probe value for all samples in the same group
aux <- data[rownames(data) %in% probe, colnames(data) %in% samples.in.group]
idx[line,3] <- mean(as.numeric(aux),na.rm = TRUE)
idx[line,4] <- median(as.numeric(aux),na.rm = TRUE)
}
# Step 2 replace
for(line in 1:nrow(idx)){
row <- idx[line,1]
col <- idx[line,2]
if(by == "mean"){
data[idx[line,1],idx[line,2]] <- idx[line,3]
} else if(by == "median") {
data[idx[line,1],idx[line,2]] <- idx[line,4]
}
}
return(data)
}
Score PanCan33 Data
Use the signature that was created and stored as ‘pcbc-stemsig.p219.rda’ to now score PanCan33 data. First load the data data.pan Replace empty values with the median probe values.
## Uses the signature to score PanCan33 data
# load TCGA 450K data (subset of 219 probes of interest)
load("data.pan.Rda")
data.pan
load("type.info.Rda") #(contains tumor type info)
type.info
testset <- replace.NA(data.pan, type.info, by="median")
======= Status Number of NA in probes ========
--------------------- Summary------------------
Min. 1st Qu. Median Mean 3rd Qu. Max.
1.00 1.00 1.00 9.36 3.00 222.00
----------- Probes with more nb of NAs -----------
cg04876127 cg10665848 cg20290850 cg01095506 cg05056545 cg15701612
222 101 35 15 11 10
===============================================
Load the signature w which should be a numeric row vector the size of the probes of interest.
load("pcbc-stemsig.p219.Rda")
w <- mm$w
w[1:10]
cg02927655 cg15948871 cg17676824 cg25552705 cg21434327 cg06678662 cg20300315 cg12257394 cg10381520 cg04193909
-0.03251136 -0.03270579 -0.03151675 -0.04281638 -0.02921891 -0.03178537 -0.02974974 -0.05360807 -0.04052402 -0.04616716
X <- testset[as.character(names(w)),]
head(X)
Convert X into a matrix
X <- as.matrix(X)
X[1:3,1:3]
TCGA-05-4384-01A-01D-1756-05 TCGA-05-4390-01A-02D-1756-05 TCGA-05-4396-01A-21D-1856-05
cg02927655 0.0978782 0.094127 0.0970064
cg15948871 0.7006290 0.675834 0.7751590
cg17676824 0.8089120 0.815386 0.8020920
Score via linear model. The resulting variable ss will be row vector 9627 in length
ss <- t(w) %*% X
ss[1,1:5]
TCGA-05-4384-01A-01D-1756-05 TCGA-05-4390-01A-02D-1756-05 TCGA-05-4396-01A-21D-1856-05 TCGA-05-4405-01A-21D-1856-05
-4.813315 -5.012777 -4.969649 -5.112913
TCGA-05-4410-01A-21D-1856-05
-4.933056
Scale the svores into a ratio from 0 to 1. and store as data frame.
## Scale the scores to be between 0 and 1
ss <- ss - min(ss)
ss <- ss / max(ss)
ss <- as.data.frame(t(ss))
colnames(ss) <- "mDNAsi"
head(ss)
Save scores to a Rda file.
save(ss, file = "TCGA_mDNAsi.Rda")
Continuous Code
# One class model - 450K data
# Files located @ https://drive.google.com/drive/folders/0BybEcxkBv6VqcmJ0QlA3bVJhc0E?usp=sharing
# Load required libraries
library(gelnet)
library(dplyr)
library(gdata)
# load PCBC metadata (contains group info)
load("pcbc.pd.f.Rda")
# load PCBC 450K data (subset of 219 probes of interest)
load("pcbc.data.Rda")
## Mean-center the data
m <- apply(pcbc.data, 1, mean )
pcbc.data.2 <- pcbc.data - m
# Define PCBC groups (SC and non.SC)
M1_smp <- pcbc.pd.f[pcbc.pd.f$Diffname_short %in% "SC",] #SC
M2_smp <- pcbc.pd.f[!(pcbc.pd.f$Diffname_short %in% "SC"),] #non-SC
# Select PCBC data
X.tr <- pcbc.data.2[, as.character(M1_smp$UID)] # 44 samples
X.bk <- pcbc.data.2[, as.character(M2_smp$UID)] # 55 samples
## Train a one-class model
mm <- gelnet(t(X.tr), NULL, 0, 1) # NULL for a one-class task
## Store the signature to a file
save( mm, file = "pcbc-stemsig.p219.Rda")
# Cross-validation with linear model:
## Perform leave-one-out cross-validation
auc <- c()
for( i in 1:ncol(X.tr) )
{
## Train a model on non-left-out data
X1 <- X.tr[,-i]
X1 <- as.matrix(X1)
K <- t(X1) %*% X1 / nrow(X1)
m1 <- gelnet.ker(K, NULL, lambda = 1)
w1 <- X1 %*% m1$v
## Score the left-out sample against the background
X.bk <- X.bk[rownames(X.tr),]
X.bk <- as.matrix(X.bk)
s.bk <- t(w1) %*% X.bk
s.bk <- unmatrix(s.bk)
s1 <- t(w1) %*% X.tr[,i]
s1 <- unmatrix(s1)
## AUC = P( left-out sample is scored above the background )
auc[i] <- sum(s1 > s.bk) / length(s.bk)
cat( "Current AUC: ", auc[i], "\n" )
cat( "Average AUC: ", mean(auc), "\n" )
}
## Uses the signature to score PanCan33 data
# load TCGA 450K data (subset of 219 probes of interest)
load("data.pan.Rda")
# Function to replace NA values with median of probe values by tumor type
source("replaceNA.R")
load("type.info.Rda") #(contains tumor type info)
testset <- replace.NA(data.pan, type.info, by = "median")
## Load the signature
load("pcbc-stemsig.p219.Rda")
w <- mm$w
X <- testset[as.character(names(w)),]
X <- as.matrix(X)
## Score via linear model
ss <- t(w) %*% X
## Scale the scores to be between 0 and 1
ss <- ss - min(ss)
ss <- ss / max(ss)
ss <- as.data.frame(t(ss))
colnames(ss) <- "mDNAsi"
save(ss, file = "TCGA_mDNAsi.Rda")
References
Sokolov, Artem, Daniel E Carlin, Evan O Paull, Robert Baertsch, and Joshua M Stuart. 2016. “Pathway-Based Genomics Prediction Using Generalized Elastic Net.” PLoS Comput Biol 12 (3). Public Library of Science: e1004790.
LS0tCnRpdGxlOiAiT25lIGNsYXNzIG1vZGVsIC0gNDUwSyBkYXRhIgphdXRob3I6ICJUYXRoaWFuZSBNYWlzdHJvIE1hbHRhIgpvdXRwdXQ6CiAgaHRtbF9ub3RlYm9vazoKICAgIGhpZ2hsaWdodDogdGFuZ28KICAgIHRoZW1lOiBqb3VybmFsCiAgICB0b2M6IHllcwogICAgdG9jX2Zsb2F0OgogICAgICBjb2xsYXBzZWQ6IG5vCiAgaHRtbF9kb2N1bWVudDoKICAgIHRvYzogeWVzCmVkaXRvcl9vcHRpb25zOgogIGNodW5rX291dHB1dF90eXBlOiBpbmxpbmUKYmlibGlvZ3JhcGh5OiBiaWJsaW9ncmFwaHkuYmliCi0tLQoqKioKIyBQdXJwb3NlIAoKKioqCgojIEludHJvZHVjdGlvbiAKCiMgU2V0dXAgCgpCZWdpbiBieSBkb3dubG9hZGluZyBhbmQgb3BlbmluZyB0aGUgZm9sbG93aW5nIHBhY2thZ2VzIGludG8geW91ciBsaWJyYXJ5LgpUaGUgW2dlbG5ldF0oaHR0cHM6Ly9jcmFuLnItcHJvamVjdC5vcmcvd2ViL3BhY2thZ2VzL2dlbG5ldC9pbmRleC5odG1sKSBsaWJyYXJ5IHdpbGwgYmUgdXNlZCB0byB0cmFpbiB0aGUgbW9kZWwsIApbZHBseXJdKGh0dHBzOi8vY3Jhbi5yLXByb2plY3Qub3JnL3dlYi9wYWNrYWdlcy9kcGx5ci9pbmRleC5odG1sKSBhbmQgW2dkYXRhXShodHRwczovL2NyYW4uci1wcm9qZWN0Lm9yZy93ZWIvcGFja2FnZXMvZ2RhdGEvaW5kZXguaHRtbCkgYXJlIHVzZWQgZm9yIGRhdGEgbWFuaXB1bGF0aW9uLgpgYGB7cn0KbGlicmFyeShnZWxuZXQpCmxpYnJhcnkoZHBseXIpCmxpYnJhcnkoZ2RhdGEpCmBgYAoKIyBEYXRhIAoKVGhlIHRoZSBwY2JjIGRhdGEgaXMgb3JnYW5pemVkIGluIGEgZGF0YSBmcmFtZSB3aXRoIDk5IHNhbXBsZXMgYXMgY29sdW1ucyBhbmQgMjE5IG1ldGh5bGF0aW9uIHByb2JlcyBhcyByb3dzLgpgYGB7cn0KbG9hZCgicGNiYy5kYXRhLlJkYSIpCnBjYmMuZGF0YQpgYGAKCmBgYHtyfQpsb2FkKCJwY2JjLnBkLmYuUmRhIikKaGVhZChwY2JjLnBkLmYpCmBgYAoKClRoZW4gZmluZCB0aGUgbWVhbiBjZW50ZXIgYnkgc3VidHJhY3RpbmcgdGhlIG1lYW4gb2YgZWFjaCBwcm9iZSBmcm9tIHRoZSBlbnRpcmUgcGNiYyBkYXRhLiAKVGhlIG1lYW4gb2YgZWFjaCBwcm9iZSBqdXN0IGJlIGluIGEgbnVtZXJpYyB2ZWN0b3IgdGhlIHNhbWUgc2l6ZSBhcyB0aGUgbnVtYmVyIG9mIHByb2JlcywgaW4gdGhpcyBjYXNlIDIxOS4KYGBge3J9Cm0gPC0gYXBwbHkocGNiYy5kYXRhLCAxLCBtZWFuKQptWzE6NV0KYGBgCmBgYHtyfQpwY2JjLmRhdGEuMiA8LSBwY2JjLmRhdGEgLSBtCnBjYmMuZGF0YS4yCmBgYAoKCklkZW50aWZ5IHN0ZW0gY2VsbHMgYW5kIGJyZWFrIHVwIGFsbCBzYW1wbGVzIGludG8gMiBncm91cHM6CgoqIFN0ZW0gY2VsbCAoJ1gudHInIG9iamVjdCkKKiBub3Qgc3RlbSBjZWxsICgnWC5iaycgb2JqZWN0KQoKYGBge3J9CiMgRGVmaW5lIFBDQkMgZ3JvdXBzIChTQyBhbmQgbm9uLlNDKQpNMV9zbXAgPC0gcGNiYy5wZC5mW3BjYmMucGQuZiREaWZmbmFtZV9zaG9ydCAlaW4lICJTQyIsXSAjU0MKTTJfc21wIDwtIHBjYmMucGQuZlshKHBjYmMucGQuZiREaWZmbmFtZV9zaG9ydCAlaW4lICJTQyIpLF0gI25vbi1TQwoKIyBTZWxlY3QgUENCQyBkYXRhClgudHIgPC0gcGNiYy5kYXRhLjJbLCBhcy5jaGFyYWN0ZXIoTTFfc21wJFVJRCldICMgNDQgc2FtcGxlcwpYLmJrIDwtIHBjYmMuZGF0YS4yWywgYXMuY2hhcmFjdGVyKE0yX3NtcCRVSUQpXSAjIDU1IHNhbXBsZXMKYGBgCgpOb3cgd2UgY2FuIGJlZ2luIHRvIHRyYWluIHRoZSB0aGUgb25lLWNsYXNzIG1vZGVsIHdpdGggdGhlIFtnZWxuZXRdKGh0dHBzOi8vd3d3LnJkb2N1bWVudGF0aW9uLm9yZy9wYWNrYWdlcy9nZWxuZXQvdmVyc2lvbnMvMS4yLjEvdG9waWNzL2dlbG5ldCkgZnVuY3Rpb24uIApUaGUgZ2VsbmV0IGZ1bmN0aW9uIGNhbiBiZSB1c2VkIGZvciBMaW5lYXIgUmVncmVzc2lvbiwgCkJpbmFyeSBDbGFzc2lmaWNhdGlvbiBhbmQgT25lIGNsYXNzIFByb2JsZW1zIGJ5IHVzaW5nIGFuIGl0ZXJhdGl2ZSBtZXRob2QgY2FsbGVkIApjb29yZGluYXRlZCBkZXNjZW50IFtAc29rb2xvdjIwMTZwYXRod2F5XS4KCmBgYHtyLCBldmFsPUZBTFNFfQpnZWxuZXQoWCwgeSwgbDEsIGwyKQpgYGAKCkl0IGhhcyBmb3VyIG1haW4gYXJndW1lbnRzIGRlc2NyaWJlZCBiZWxvdzoKCiogKipYKio6IG4gYnkgcCBtYXRyaXggPT4gdHJhbnNwb3NlKCdYLnInKQoqICoqeSoqOiAnTlVMTCcgZm9yIG9uZSBjbGFzcyBtb2RlbHMgCiogKipsMSoqOiBjb2VmZmljaWVudCBmb3IgdGhlIEwxLW5vcm0gcGVuYWx0eSA9PiAnMCcgCiogKipsMioqOiBjb2VmZmljaWVudCBmb3IgdGhlIEwyLW5vcm0gcGVuYWx0eSA9PiAnMScKCk1ha2Ugc3VyZSB5b3UgdHJhbnNwb3NlIHRoZSBtYXRyaXggc28gdGhhdCB0aGUgZ2VuZXMgYXJlIGxpc3RlZCBhcyByb3dzIGFuZCBzYW1wbGVzIGFzIGNvbHVtbnMuClRoZW4gc3RvcmUgdGhlIHJlc3VsdHMgYXMgYSB0c3YgZmlsZSAocGNiYy1zdGVtc2lnLnAyMTkucmRhKS4gCmBgYHtyfQojIyBUcmFpbiBhIG9uZS1jbGFzcyBtb2RlbAptbSA8LSBnZWxuZXQodChYLnRyKSwgTlVMTCwgMCwgMSkgI05VTEwgZm9yIGEgb25lLWNsYXNzIHRhc2sgCgojIyBTdG9yZSB0aGUgc2lnbmF0dXJlIHRvIGEgZmlsZQpzYXZlKCBtbSwgZmlsZSA9ICJwY2JjLXN0ZW1zaWcucDIxOS5SZGEiKQpgYGAKCgoKIyBMZWF2ZSBPbmUgT3V0IENyb3NzLVZhbGlkYXRpb24gCgpUbyB0ZXN0IGhvdyB0aGUgbW9kZWwncyBwZXJmb3JtYW5jZSBieSB1c2luZyBsZWF2ZSBvbmUgb3V0IGNyb3NzLXZhbGlkYXRpb24uIApUaGlzIHByb2Nlc3MgaGFzIHRocmVlIHN0ZXBzOiAKCjEuIFRyYWluIG1vZGVsIG9uIG5vbi1sZWZ0LW91dCBkYXRhCjIuIFNjb3JlIHRoZSBsZWZ0LW91dCBzYW1wbGUgYWdhaW5zdCB0aGUgYmFja2dyb3VuZAozLiBBVUMgPSBQKCBsZWZ0LW91dCBzYW1wbGUgaXMgc2NvcmVkIGFib3ZlIHRoZSBiYWNrZ3JvdW5kICkKCmBgYHtyLCByZXN1bHRzPSJoaWRlIn0KIyBDcm9zcy12YWxpZGF0aW9uIHdpdGggbGluZWFyIG1vZGVsOgojIFBlcmZvcm0gbGVhdmUtb25lLW91dCBjcm9zcy12YWxpZGF0aW9uCmF1YyA8LSBjKCkKZm9yKGkgaW4gMTpuY29sKFgudHIpKSB7CiAgIyMgVHJhaW4gYSBtb2RlbCBvbiBub24tbGVmdC1vdXQgZGF0YQogIFgxIDwtIFgudHJbLC1pXQogIFgxIDwtIGFzLm1hdHJpeChYMSkKICBLIDwtIHQoWDEpICUqJSBYMSAvIG5yb3coWDEpCiAgbTEgPC0gZ2VsbmV0LmtlcihLLCBOVUxMLCBsYW1iZGEgPSAxKQogIHcxIDwtIFgxICUqJSBtMSR2CiAgCiAgIyMgU2NvcmUgdGhlIGxlZnQtb3V0IHNhbXBsZSBhZ2FpbnN0IHRoZSBiYWNrZ3JvdW5kCiAgWC5iayA8LSBYLmJrW3Jvd25hbWVzKFgudHIpLF0KICBYLmJrIDwtIGFzLm1hdHJpeChYLmJrKQogIHMuYmsgPC0gdCh3MSkgJSolIFguYmsKICBzLmJrIDwtIHVubWF0cml4KHMuYmspCiAgCiAgczEgPC0gdCh3MSkgJSolIFgudHJbLGldCiAgczEgPC0gdW5tYXRyaXgoczEpCiAgCiAgIyMgQVVDID0gUCggbGVmdC1vdXQgc2FtcGxlIGlzIHNjb3JlZCBhYm92ZSB0aGUgYmFja2dyb3VuZCApCiAgYXVjW2ldIDwtIHN1bShzMSA+IHMuYmspIC8gbGVuZ3RoKHMuYmspCiAgY2F0KCAiQ3VycmVudCBBVUM6ICIsIGF1Y1tpXSwgIlxuIiApCiAgY2F0KCAiQXZlcmFnZSBBVUM6ICIsIG1lYW4oYXVjKSwgIlxuIiApCn0KYGBgCgoKSWYgdGhlIHZhbGlkYXRpb24gaXMgc3VjY2Vzc2Z1bCB5b3Ugd2lsbCBub3RpY2UgdGhhdCB0aGUgYXVjIHZhcmlhYmxlIHdpbGwgYmUgYSBudW1lcmljIHZlY3RvciBjb25zaXN0aW5nIG9mIG9ubHkgMSdzLiAgCmBgYHtyfQpwcmludChzMSkKYGBgCmBgYHtyfQpoZWFkKGF1YykKYWxsKGF1YyA9PSAxKQpgYGAKCiMgUmVwbGFjZSBOQSAKClRoZSBSZXBsYWNlIE5BIGZ1bmN0aW9uIGlzIHVzZWQgdG8gcmVwbGFjZSBhbnkgdmFsdWVzIHRoYXQgYXJlIE5BIChub3QgYXZhaWxhYmxlKSAKd2l0aCAgZWl0aGVyIHRoZSBtZWFuIG9yIHRoZSBtZWRpYW4gdmFsdWUgb2YgdGhlIHByb2JlIGZvciBhIGdpdmUgZ3JvdXAuIAoKMS4gY2hlY2sgZm9yIE5BIHZhbHVlcwoyLiBsb2NhdGUgdGhlIE5BIHZhbHVlcyAKMy4gY2FsY3VsYXRlIHRoZSBtZWFuIG9yIG1lZGlhbiBmb3IgdGhlIHByb2JlIHdoZXJlIGVhY2ggTkEgaXMgZm91bmQgCjQuIHJlcGxhY2UgdGhlIHZhbHVlcwoKYGBge3IsIG1lc3NhZ2UgPSBGQUxTRSwgZXJyb3IgPSBGQUxTRX0KcmVwbGFjZS5OQSA8LWZ1bmN0aW9uKGRhdGEsdHlwZS5pbmZvLGJ5ID0gIm1lYW4iKXsKICBpZighImdyb3VwIiAlaW4lIGNvbG5hbWVzKHR5cGUuaW5mbykpIHN0b3AoInR5cGUuaW5mbyBtdXN0IGhhdmUgZ3JvdXAgY29sdW1uIikKICBpZighInNhbXBsZSIgJWluJSBjb2xuYW1lcyh0eXBlLmluZm8pKSBzdG9wKCJ0eXBlLmluZm8gbXVzdCBoYXZlIGEgc2FtcGxlIGNvbHVtbiIpCiAgCiAgIyBEbyB3ZSBoYXZlIE5Bcz8KICBpZihpcy5uYSh0YWJsZShpcy5uYShkYXRhKSlbIlRSVUUiXSkpewogICAgbWVzc2FnZSgiTm8gTkFzIHdlcmUgZm91bmQiKQogICAgcmV0dXJuKGRhdGEpCiAgfQogICMgZ2V0IE5BcyBpbmRleCAKICBpZHggPC0gd2hpY2goaXMubmEoZGF0YSkgPT0gVFJVRSxhcnIuaW5kPVRSVUUpCiAgY291bnQgPC0gdGFibGUocm93bmFtZXMoaWR4KSkKICBtZXNzYWdlKCI9PT09PT09IFN0YXR1cyBOdW1iZXIgb2YgTkEgaW4gcHJvYmVzID09PT09PT09IikKICBtZXNzYWdlKCItLS0tLS0tLS0tLS0tLS0tLS0tLS0gU3VtbWFyeS0tLS0tLS0tLS0tLS0tLS0tLSIpCiAgcHJpbnQoc3VtbWFyeShhcy5udW1lcmljKGNvdW50KSkpCiAgbWVzc2FnZSgiXG4tLS0tLS0tLS0tLSBQcm9iZXMgd2l0aCBtb3JlIG5iIG9mIE5BcyAtLS0tLS0tLS0tLSIpCiAgcHJpbnQoaGVhZChzb3J0KGNvdW50LGRlY3JlYXNpbmcgPSBUKSkpCiAgbWVzc2FnZSgiPT09PT09PT09PT09PT09PT09PT09PT09PT09PT09PT09PT09PT09PT09PT09PT0iKQogICAgICAgICAgICAgICAgIAogIGlkeCA8LSBjYmluZChpZHgsIG1lYW4gPSBOQSwgbWVkaWFuID0gTkEpCiAgCiAgIyBGb3IgZWFjaCBOQSB2YWx1ZSBjYWxjdWxhdGUgdGhlIG1lYW4gZm9yIHRoZSBzYW1lIHByb2JlIGZvciB0aGUgc2FtcGxlcwogICMgd2hlcmUgaXQgYmVsb25ncwogIGZvcihsaW5lIGluIDE6bnJvdyhpZHgpKXsKICAgIHJvdyA8LSBpZHhbbGluZSwxXQogICAgY29sIDwtIGlkeFtsaW5lLDJdCiAgICBwcm9iZSA8LSByb3duYW1lcyhpZHgpW2xpbmVdCiAgICBzYW1wbGUgPC0gY29sbmFtZXMoZGF0YSlbY29sXQogICAgZ3JvdXAgPC0gdHlwZS5pbmZvW3R5cGUuaW5mbyRzYW1wbGUgPT0gc2FtcGxlLCJncm91cCJdCiAgICBzYW1wbGVzLmluLmdyb3VwIDwtIHR5cGUuaW5mb1t0eXBlLmluZm8kZ3JvdXAgPT0gZ3JvdXAsXSRzYW1wbGUKICAgIAogICAgIyBnZXQgdGhlIHByb2JlIHZhbHVlIGZvciBhbGwgc2FtcGxlcyBpbiB0aGUgc2FtZSBncm91cCAKICAgIGF1eCA8LSBkYXRhW3Jvd25hbWVzKGRhdGEpICVpbiUgcHJvYmUsIGNvbG5hbWVzKGRhdGEpICVpbiUgc2FtcGxlcy5pbi5ncm91cF0gCiAgICAKICAgIGlkeFtsaW5lLDNdIDwtIG1lYW4oYXMubnVtZXJpYyhhdXgpLG5hLnJtID0gVFJVRSkKICAgIGlkeFtsaW5lLDRdIDwtIG1lZGlhbihhcy5udW1lcmljKGF1eCksbmEucm0gPSBUUlVFKQogIH0KICAjIFN0ZXAgMiByZXBsYWNlCiAgZm9yKGxpbmUgaW4gMTpucm93KGlkeCkpewogICAgcm93IDwtIGlkeFtsaW5lLDFdCiAgICBjb2wgPC0gaWR4W2xpbmUsMl0KICAgIGlmKGJ5ID09ICJtZWFuIil7CiAgICAgIGRhdGFbaWR4W2xpbmUsMV0saWR4W2xpbmUsMl1dIDwtIGlkeFtsaW5lLDNdICAKICAgIH0gZWxzZSBpZihieSA9PSAibWVkaWFuIikgeyAKICAgICAgZGF0YVtpZHhbbGluZSwxXSxpZHhbbGluZSwyXV0gPC0gaWR4W2xpbmUsNF0KICAgIH0KICB9CiAgcmV0dXJuKGRhdGEpCn0KYGBgCgojIFNjb3JlIFBhbkNhbjMzIERhdGEgCgpVc2UgdGhlIHNpZ25hdHVyZSB0aGF0IHdhcyBjcmVhdGVkIGFuZCBzdG9yZWQgYXMgJ3BjYmMtc3RlbXNpZy5wMjE5LnJkYScgdG8gbm93IHNjb3JlIFBhbkNhbjMzIGRhdGEuIApGaXJzdCBsb2FkIHRoZSBkYXRhIGRhdGEucGFuClJlcGxhY2UgZW1wdHkgdmFsdWVzIHdpdGggdGhlIG1lZGlhbiBwcm9iZSB2YWx1ZXMuCgpgYGB7cn0KIyMgVXNlcyB0aGUgc2lnbmF0dXJlIHRvIHNjb3JlIFBhbkNhbjMzIGRhdGEKIyBsb2FkIFRDR0EgNDUwSyBkYXRhIChzdWJzZXQgb2YgMjE5IHByb2JlcyBvZiBpbnRlcmVzdCkKbG9hZCgiZGF0YS5wYW4uUmRhIikgCmRhdGEucGFuCmBgYApgYGB7cn0KbG9hZCgidHlwZS5pbmZvLlJkYSIpICMoY29udGFpbnMgdHVtb3IgdHlwZSBpbmZvKQp0eXBlLmluZm8KYGBgCmBgYHtyfQp0ZXN0c2V0IDwtIHJlcGxhY2UuTkEoZGF0YS5wYW4sIHR5cGUuaW5mbywgYnk9Im1lZGlhbiIpIApoZWFkKHRlc3RzZXQpCmBgYAoKCkxvYWQgdGhlIHNpZ25hdHVyZSB3IHdoaWNoIHNob3VsZCBiZSBhIG51bWVyaWMgcm93IHZlY3RvciB0aGUgc2l6ZSBvZiB0aGUgcHJvYmVzIG9mIGludGVyZXN0LiAKYGBge3J9CmxvYWQoInBjYmMtc3RlbXNpZy5wMjE5LlJkYSIpCncgPC0gbW0kdwp3WzE6MTBdCmBgYAoKCmBgYHtyfQpYIDwtIHRlc3RzZXRbYXMuY2hhcmFjdGVyKG5hbWVzKHcpKSxdCmhlYWQoWCkKYGBgCgpDb252ZXJ0IFggaW50byBhIG1hdHJpeApgYGB7cn0KWCA8LSBhcy5tYXRyaXgoWCkKWFsxOjMsMTozXQpgYGAKCgpTY29yZSB2aWEgbGluZWFyIG1vZGVsLiBUaGUgcmVzdWx0aW5nIHZhcmlhYmxlIHNzIHdpbGwgYmUgcm93IHZlY3RvciA5NjI3IGluIGxlbmd0aCAKYGBge3J9CnNzIDwtIHQodykgJSolIFgKc3NbMSwxOjVdCmBgYAoKU2NhbGUgdGhlIHN2b3JlcyBpbnRvIGEgcmF0aW8gZnJvbSAwIHRvIDEuIGFuZCBzdG9yZSBhcyBkYXRhIGZyYW1lLiAKYGBge3J9CiMjIFNjYWxlIHRoZSBzY29yZXMgdG8gYmUgYmV0d2VlbiAwIGFuZCAxCnNzIDwtIHNzIC0gbWluKHNzKQpzcyA8LSBzcyAvIG1heChzcykKc3MgPC0gYXMuZGF0YS5mcmFtZSh0KHNzKSkKY29sbmFtZXMoc3MpIDwtICJtRE5Bc2kiICAgCmhlYWQoc3MpCmBgYAoKU2F2ZSBzY29yZXMgdG8gYSBSZGEgZmlsZS4KYGBge3J9CnNhdmUoc3MsIGZpbGUgPSAiVENHQV9tRE5Bc2kuUmRhIikKYGBgCgojIENvbnRpbnVvdXMgQ29kZSAKCmBgYHtyLGV2YWwgPSBGQUxTRX0KIyBPbmUgY2xhc3MgbW9kZWwgLSA0NTBLIGRhdGEKIyBGaWxlcyBsb2NhdGVkIEAgaHR0cHM6Ly9kcml2ZS5nb29nbGUuY29tL2RyaXZlL2ZvbGRlcnMvMEJ5YkVjeGtCdjZWcWNtSjBRbEEzYlZKaGMwRT91c3A9c2hhcmluZwoKIyBMb2FkIHJlcXVpcmVkIGxpYnJhcmllcwpsaWJyYXJ5KGdlbG5ldCkKbGlicmFyeShkcGx5cikKbGlicmFyeShnZGF0YSkKCiMgbG9hZCBQQ0JDIG1ldGFkYXRhIChjb250YWlucyBncm91cCBpbmZvKQpsb2FkKCJwY2JjLnBkLmYuUmRhIikKIyBsb2FkIFBDQkMgNDUwSyBkYXRhIChzdWJzZXQgb2YgMjE5IHByb2JlcyBvZiBpbnRlcmVzdCkKbG9hZCgicGNiYy5kYXRhLlJkYSIpCgojIyBNZWFuLWNlbnRlciB0aGUgZGF0YQptIDwtIGFwcGx5KHBjYmMuZGF0YSwgMSwgbWVhbiApCnBjYmMuZGF0YS4yIDwtIHBjYmMuZGF0YSAtIG0KCiMgRGVmaW5lIFBDQkMgZ3JvdXBzIChTQyBhbmQgbm9uLlNDKQpNMV9zbXAgPC0gcGNiYy5wZC5mW3BjYmMucGQuZiREaWZmbmFtZV9zaG9ydCAlaW4lICJTQyIsXSAjU0MKTTJfc21wIDwtIHBjYmMucGQuZlshKHBjYmMucGQuZiREaWZmbmFtZV9zaG9ydCAlaW4lICJTQyIpLF0gI25vbi1TQwoKIyBTZWxlY3QgUENCQyBkYXRhClgudHIgPC0gcGNiYy5kYXRhLjJbLCBhcy5jaGFyYWN0ZXIoTTFfc21wJFVJRCldICMgNDQgc2FtcGxlcwpYLmJrIDwtIHBjYmMuZGF0YS4yWywgYXMuY2hhcmFjdGVyKE0yX3NtcCRVSUQpXSAjIDU1IHNhbXBsZXMKCiMjIFRyYWluIGEgb25lLWNsYXNzIG1vZGVsCm1tIDwtIGdlbG5ldCh0KFgudHIpLCBOVUxMLCAwLCAxKSAjIE5VTEwgZm9yIGEgb25lLWNsYXNzIHRhc2sgCgojIyBTdG9yZSB0aGUgc2lnbmF0dXJlIHRvIGEgZmlsZQpzYXZlKCBtbSwgZmlsZSA9ICJwY2JjLXN0ZW1zaWcucDIxOS5SZGEiKQoKIyBDcm9zcy12YWxpZGF0aW9uIHdpdGggbGluZWFyIG1vZGVsOgojIyBQZXJmb3JtIGxlYXZlLW9uZS1vdXQgY3Jvc3MtdmFsaWRhdGlvbgphdWMgPC0gYygpCmZvciggaSBpbiAxOm5jb2woWC50cikgKQp7CiAgIyMgVHJhaW4gYSBtb2RlbCBvbiBub24tbGVmdC1vdXQgZGF0YQogIFgxIDwtIFgudHJbLC1pXQogIFgxIDwtIGFzLm1hdHJpeChYMSkKICBLIDwtIHQoWDEpICUqJSBYMSAvIG5yb3coWDEpCiAgbTEgPC0gZ2VsbmV0LmtlcihLLCBOVUxMLCBsYW1iZGEgPSAxKQogIHcxIDwtIFgxICUqJSBtMSR2CiAgCiAgIyMgU2NvcmUgdGhlIGxlZnQtb3V0IHNhbXBsZSBhZ2FpbnN0IHRoZSBiYWNrZ3JvdW5kCiAgWC5iayA8LSBYLmJrW3Jvd25hbWVzKFgudHIpLF0KICBYLmJrIDwtIGFzLm1hdHJpeChYLmJrKQogIHMuYmsgPC0gdCh3MSkgJSolIFguYmsKICBzLmJrIDwtIHVubWF0cml4KHMuYmspCiAgCiAgczEgPC0gdCh3MSkgJSolIFgudHJbLGldCiAgczEgPC0gdW5tYXRyaXgoczEpCiAgCiAgIyMgQVVDID0gUCggbGVmdC1vdXQgc2FtcGxlIGlzIHNjb3JlZCBhYm92ZSB0aGUgYmFja2dyb3VuZCApCiAgYXVjW2ldIDwtIHN1bShzMSA+IHMuYmspIC8gbGVuZ3RoKHMuYmspCiAgY2F0KCAiQ3VycmVudCBBVUM6ICIsIGF1Y1tpXSwgIlxuIiApCiAgY2F0KCAiQXZlcmFnZSBBVUM6ICIsIG1lYW4oYXVjKSwgIlxuIiApCn0KCiMjIFVzZXMgdGhlIHNpZ25hdHVyZSB0byBzY29yZSBQYW5DYW4zMyBkYXRhCiMgbG9hZCBUQ0dBIDQ1MEsgZGF0YSAoc3Vic2V0IG9mIDIxOSBwcm9iZXMgb2YgaW50ZXJlc3QpCmxvYWQoImRhdGEucGFuLlJkYSIpIAojIEZ1bmN0aW9uIHRvIHJlcGxhY2UgTkEgdmFsdWVzIHdpdGggbWVkaWFuIG9mIHByb2JlIHZhbHVlcyBieSB0dW1vciB0eXBlCnNvdXJjZSgicmVwbGFjZU5BLlIiKQpsb2FkKCJ0eXBlLmluZm8uUmRhIikgIyhjb250YWlucyB0dW1vciB0eXBlIGluZm8pCnRlc3RzZXQgPC0gcmVwbGFjZS5OQShkYXRhLnBhbiwgdHlwZS5pbmZvLCBieSA9ICJtZWRpYW4iKSAKCiMjIExvYWQgdGhlIHNpZ25hdHVyZQpsb2FkKCJwY2JjLXN0ZW1zaWcucDIxOS5SZGEiKQp3IDwtIG1tJHcKClggPC0gdGVzdHNldFthcy5jaGFyYWN0ZXIobmFtZXModykpLF0KWCA8LSBhcy5tYXRyaXgoWCkKCiMjIFNjb3JlIHZpYSBsaW5lYXIgbW9kZWwKc3MgPC0gdCh3KSAlKiUgWAojIyBTY2FsZSB0aGUgc2NvcmVzIHRvIGJlIGJldHdlZW4gMCBhbmQgMQpzcyA8LSBzcyAtIG1pbihzcykKc3MgPC0gc3MgLyBtYXgoc3MpCnNzIDwtIGFzLmRhdGEuZnJhbWUodChzcykpCgpjb2xuYW1lcyhzcykgPC0gIm1ETkFzaSIgICAKc2F2ZShzcywgZmlsZSA9ICJUQ0dBX21ETkFzaS5SZGEiKQpgYGAKCiMgUmVmZXJlbmNlcwoKCgoKCgo=